Here are the details of the Bayesian multi-level modelling. Our general approach is:
Define Priors
In this section we will specify some priors. We then then use a prior-predictive check to assess whether the prior is reasonable or not (i.e., on the same order of magnitude as our measurements).
Fixed Effects
Our independent variable is a categorical factor with 16 levels. We will drop the intercept from our model and instead fit a coefficent for each factor level (\(y \sim x - 0\)). As our dependant variable will be log-transformed, we can use the priors below:
prior <- c(
set_prior("normal(0,2)", class = "b"),
set_prior("cauchy(0,2)", class = "sigma"))
Group-level Effects
We will keep the default weakly informative priors for the group-level (‘random’) effects. From the brms documentation:
[…] restricted to be non-negative and, by default, have a half student-t prior with 3 degrees of freedom and a scale parameter that depends on the standard deviation of the response after applying the link function. Minimally, the scale parameter is 10. This prior is used (a) to be only very weakly informative in order to influence results as few as possible, while (b) providing at least some regularization to considerably improve convergence and sampling efficiency.
Prior Predictive Check
Now we can specify our Bayesian multi-level model and priors. Note that as we are using sample_prior = 'only', the model will not learn anything from our data.
m_prior <- brm(data = d_eeg,
rms ~ wg-1 + (1|subject),
family = "lognormal",
prior = prior,
iter = n_iter ,
sample_prior = 'only')
We can use this model to generate data.
## Observations: 320,000
## Variables: 2
## $ key <chr> "P2", "P2", "P2", "P2", "P2", "P2", "P2", "P2", "P2", "P2"…
## $ value <dbl> 0.11096794, 0.99479785, 0.48693538, 0.03488185, 7.09463726…
We can see that i) our priors are relatively weak as the predictions span several orders of magnitide, and ii) our empirical data falls within this range.
Compute Posterior
Fit Model to EEG Data
We will now fit the model to the data.
## Family: lognormal
## Links: mu = identity; sigma = identity
## Formula: rms ~ wg - 1 + (1 | subject)
## Data: d_eeg (Number of observations: 400)
## Samples: 4 chains, each with iter = 10000; warmup = 5000; thin = 1;
## total post-warmup samples = 20000
##
## Group-Level Effects:
## ~subject (Number of levels: 25)
## Estimate Est.Error l-95% CI u-95% CI Rhat Bulk_ESS Tail_ESS
## sd(Intercept) 0.38 0.06 0.28 0.52 1.01 1151 2036
##
## Population-Level Effects:
## Estimate Est.Error l-95% CI u-95% CI Rhat Bulk_ESS Tail_ESS
## wgCM -0.49 0.09 -0.66 -0.31 1.01 803 1575
## wgCMM -0.12 0.09 -0.30 0.05 1.01 834 1565
## wgP2 -0.95 0.09 -1.13 -0.78 1.01 792 1539
## wgP3 -0.82 0.09 -1.00 -0.65 1.01 778 1372
## wgP31M -0.43 0.09 -0.61 -0.26 1.00 793 1280
## wgP3M1 -0.14 0.09 -0.32 0.03 1.01 798 1494
## wgP4 -0.48 0.09 -0.66 -0.31 1.00 817 1387
## wgP4G -0.26 0.09 -0.43 -0.08 1.01 800 1515
## wgP4M 0.15 0.09 -0.03 0.33 1.01 823 1490
## wgP6 -0.48 0.09 -0.67 -0.31 1.01 802 1564
## wgP6M -0.05 0.09 -0.22 0.13 1.01 818 1413
## wgPG -0.93 0.09 -1.11 -0.76 1.01 808 1422
## wgPGG -0.73 0.09 -0.90 -0.55 1.01 788 1543
## wgPM -0.59 0.09 -0.77 -0.42 1.00 816 1459
## wgPMG -0.26 0.09 -0.44 -0.09 1.01 798 1471
## wgPMM -0.07 0.09 -0.25 0.10 1.01 791 1532
##
## Family Specific Parameters:
## Estimate Est.Error l-95% CI u-95% CI Rhat Bulk_ESS Tail_ESS
## sigma 0.23 0.01 0.21 0.25 1.00 8335 10925
##
## Samples were drawn using sampling(NUTS). For each parameter, Bulk_ESS
## and Tail_ESS are effective sample size measures, and Rhat is the potential
## scale reduction factor on split chains (at convergence, Rhat = 1).
We will now look at the model’s predicts for the average participant (i.e, ignoring the random intercepts).
Fit Model to Psychophysical Data
We will now fit the model to the data.
## Family: lognormal
## Links: mu = identity; sigma = identity
## Formula: threshold ~ wg - 1 + (1 | subject)
## Data: d_dispthresh (Number of observations: 186)
## Samples: 4 chains, each with iter = 10000; warmup = 5000; thin = 1;
## total post-warmup samples = 20000
##
## Group-Level Effects:
## ~subject (Number of levels: 12)
## Estimate Est.Error l-95% CI u-95% CI Rhat Bulk_ESS Tail_ESS
## sd(Intercept) 0.39 0.11 0.23 0.65 1.00 3014 5017
##
## Population-Level Effects:
## Estimate Est.Error l-95% CI u-95% CI Rhat Bulk_ESS Tail_ESS
## wgCM -0.22 0.17 -0.54 0.12 1.00 2375 3804
## wgCMM -1.15 0.17 -1.48 -0.81 1.00 2304 3857
## wgP2 -0.29 0.17 -0.62 0.06 1.00 2404 4013
## wgP3 -0.69 0.17 -1.00 -0.35 1.00 2319 4007
## wgP31M -1.18 0.17 -1.50 -0.84 1.00 2404 4055
## wgP3M1 -1.35 0.17 -1.67 -1.00 1.00 2348 4000
## wgP4 -0.89 0.17 -1.21 -0.55 1.00 2401 4169
## wgP4G -1.22 0.17 -1.55 -0.87 1.00 2377 4288
## wgP4M -1.29 0.17 -1.61 -0.95 1.00 2320 3621
## wgP6 -1.20 0.17 -1.53 -0.86 1.00 2417 4441
## wgP6M -1.42 0.17 -1.74 -1.07 1.00 2395 3597
## wgPG 0.36 0.17 0.03 0.71 1.00 2412 4070
## wgPGG -0.31 0.17 -0.64 0.03 1.00 2300 4383
## wgPM -0.79 0.17 -1.11 -0.44 1.00 2340 3847
## wgPMG -0.97 0.17 -1.29 -0.62 1.00 2369 3987
## wgPMM -1.17 0.17 -1.49 -0.83 1.00 2270 4049
##
## Family Specific Parameters:
## Estimate Est.Error l-95% CI u-95% CI Rhat Bulk_ESS Tail_ESS
## sigma 0.42 0.02 0.37 0.46 1.00 12058 12526
##
## Samples were drawn using sampling(NUTS). For each parameter, Bulk_ESS
## and Tail_ESS are effective sample size measures, and Rhat is the potential
## scale reduction factor on split chains (at convergence, Rhat = 1).
EEG Control Data
We will also fit models to the control data. As we can see from Figure 2.4, the group differences are much smaller.
## Warning: Bulk Effective Samples Size (ESS) is too low, indicating posterior means and medians may be unreliable.
## Running the chains for more iterations may help. See
## http://mc-stan.org/misc/warnings.html#bulk-ess
